Integrated transcriptome and endogenous hormone analysis provides new insights into callus proliferation in Osmanthus fragrans

Osmanthus fragrans is an important evergreen species with both medicinal and ornamental value in China. Given the low efficiency of callus proliferation and the difficulty of adventitious bud differentiation, tissue culture and regeneration systems have not been successfully established for this species. To understand the mechanism of callus proliferation, transcriptome sequencing and endogenous hormone content determination were performed from the initial growth stages to the early stages of senescence on O. fragrans calli. In total, 47,340 genes were identified by transcriptome sequencing, including 1798 previously unidentified genes specifically involved in callus development. Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis of differentially expressed genes (DEGs) was significantly enriched in plant hormone signal transduction pathways. Furthermore, our results from the orthogonal projections to latent structures discrimination analysis (OPLS-DA) of six typical hormones in five development stages of O. fragrans calli showed jasmonic acid (JA) could play important role in the initial stages of calli growth, whereas JA and auxin (IAA) were dominant in the early stages of calli senescence. Based on the weighted gene co-expression network analysis, OfSRC2, OfPP2CD5 and OfARR1, OfPYL3, OfEIL3b were selected as hub genes from the modules with the significant relevance to JA and IAA respectively. The gene regulation network and quantitative real-time PCR implied that during the initial stages of callus growth, the transcription factors (TFs) OfERF4 and OfMYC2a could down-regulate the expression of hub genes OfSRC2 and OfPP2CD5, resulting in decreased JA content and rapid callus growth; during the late stage of callus growth, the TFs OfERF4, OfMYC2a and OfTGA21c, OfHSFA1 could positively regulate the expression of hub genes OfSRC2, OfPP2CD5 and OfARR1, OfPYL3, OfEIL3b, respectively, leading to increased JA and IAA contents and inducing the senescence of O. fragrans calli. Hopefully, our results could provide new insights into the molecular mechanism of the proliferation of O. fragrans calli.

IAA, GA 3 , ABA, BR, and ZR, were determined in the callus of O. fragrans at different development stages by enzyme linked immunosorbent assay (ELISA) (Fig. 1A). With increases in growing days, the GA 3 content in the callus of O. fragrans first decreased and then increased, with the lowest content at 45 days and the highest at 75 days. The IAA content in calli at 25 days, 45 days, 55 days, and 65 days was unchanged, although the content suddenly increased at 75 days, which was significantly higher than in the other stages. The BR content peaked at 65 days of growth, and then decreased, being significantly different from that at 25 days, 45 days, and 55 days, but not at 75 days. The ABA content in calli was high at all five developmental stages examined, peaking at 121.12 ng/g·FW at 25 days, which was significantly different from that at other stages. The expression of JA and ZR first decreased, then increased, and then decreased again. The JA content had increased by 55 days, decreased at 65 days and increased again at 75 days, whereas the ZR content was decreased at 55 days, increased at 65 days, and decreased at 75 days.
Simca14.1 software was used for Orthogonal projections to latent structures discrimination analysis (OPLS-DA) of the six endogenous hormones measured in the callus of O. fragrans. Data points for the six hormones fell within the 95% confidence interval (CI), and the biological repeated data points in the same period had relatively clear clustering, with data points associated with the different developmental stages clustering in different regions (Fig. 1B). According to the highest value of variable important in projection (VIP) , JA was the main endogenous hormone during the initial stages of callus growth, whereas JA and IAA were the major endogenous hormones during the early stages of callus senescence (Table S1).
Transcriptome sequencing results. In  www.nature.com/scientificreports/ more. In this study, 1798 previously unidentified genes were obtained by comparing our data with the whole genome of O. fragrans. Differentially expressed genes (DEGs) during the callus growth were compared, revealing 3164 significantly different genes between 25 days versus 45 days, including 898 upregulated genes and 2266 downregulated genes (Fig. 2). However, there were 1895 significantly differently expressed genes between 65 days versus 75 days, of which 543 were upregulated and 1316 were downregulated (Fig. 2). www.nature.com/scientificreports/ biosynthesis of phenylpropanol ( Fig. S3). At 65 days versus 75 days, 396 genes were enriched in 100 metabolic pathways, with significant enrichment in eight metabolic pathways, including plant hormone signaling transduction, phenylpropanol biosynthesis, and plant circadian rhythm (Fig. S4).

Identification and functional annotation of
Identification of DEG co-expression modules by WGCNA. RNA-seq sequencing was used to analyze the DEGs of callus proliferation. To further identify the hub genes that regulate the proliferation of O. fragrans calli, WGCNA was also performed on the endogenous content hormones and the transcriptome (Fig. 3). The IDs of selected genes were shown in www.nature.com/scientificreports/ PCR validation. The trends in qRT-PCR expression were consistent with those of the FPKM values ( Fig. 4A). However, some genes showed different trends, possibly because of differences between expression detection in qRT-PCR versus transcriptome sequencing. Further correlation analysis showed that a correlation coefficient of R 2 = 0.6321, indicating that the transcriptome data were accurate and reliable, and could be used in subsequent experiments (Fig. 4B).

The co-expression network reveals the regulatory network of endogenous hormones JA and IAA. To understand the regulatory network relationship between hub genes and transcription factors (TFs)
in JA and IAA, TFs corresponding to two key genes in 'royalblue' , the upregulated module corresponding to JA, were screened. In addition, in 'brown' , the upregulated module corresponding to IAA, TFs corresponding to six hub genes were screened, and network regulation maps were drawn. Visualized in Cytoscape, nine nodes in the endogenous hormone JA regulatory network were connected to nine edges (Fig. 5A), and 35 nodes in the IAA regulatory network were connected to 55 edges (Fig. 5B). Then, 12 genes, including hub genes and corresponding TFs, were analyzed by qRT-PCR ( Fig. 6A) using the primers listed in Table S3. Using R 2 > 0.6 as the standard, five pairs of combinations were screened out with interaction relationships: OfMYC2a and OfPP2CD5, OfERF4 and OfSRC2, OfHSFA1 and OfEIL3b, OfTGA21c and OfARR1, and OfHSFA1 and OfPYL3 (Fig. 6B). Given these results, we propose a model based on changes in hormone content and network regulation maps during callus proliferation in O. fragrans (Fig. 7). During the initial growth stages, TFs (OfERF4 and OfMYC2a) could negatively regulate the hub genes (OfSRC2 and OfPP2CD5), decreasing the JA content and leading to rapid callus growth. During the early stages of calli senescence, the TFs OfTGA21c and OfHSFA1 www.nature.com/scientificreports/ could positively regulate hub genes (OfARR1 and OfPYL3, OfEIL3b) and the hub genes (OfSRC2 and OfPP2CD5) could also regulated by TFs OfERF4 and OfMYC2a respectively, leading to increased IAA and JA contents and the senescence of O. fragrans calli.

Discussion
The amount and type of endogenous hormones have an important influence on the growth and development of plants 27 , such as seed germination 28 , cell differentiation 29 , root growth 30 , and fruit ripening and shedding 31 . In the tissue culture of O. fragrans calli, we find that the calli experienced from the initial growth stages to the early stages of calli senescence. JA and IAA are important hormones involved in the regulation of many physiological processes of plant growth and development [32][33][34] . The five growth stages of callus formation were divided according to the formation time and development status of O. fragrans callus (Fig. 8). Each growth stage was regulated by a variety of endogenous hormones. The content of endogenous hormones in calli of O. fragrans during these stages www.nature.com/scientificreports/ was determined. During the growth period from 25 to 45 days, the JA content decreased gradually (Fig. 1A), and the calli grew rapidly at this time, indicating that a low concentration of JA promotes the rapid proliferation of O. fragrans calli, similar to results of studies on callus proliferation in Allium sativum 35 . In general, IAA can promote the growth of plant callus. In this study, the content of IAA remained unchanged from 25 to 65 days (Fig. 1A), but increased suddenly from 65 to 75 days, suggesting that the increased IAA concentration inhibited callus growth. The possible reason was that a certain level of IAA could promote the proliferation of callus, while higher levels of IAA could exert an inhibitory effect on this process, which was consistent with the results of previous studies 36 . These results suggested that the changes of JA and IAA contents affected the calli growth of O. fragrans. Endogenous hormones can coregulate the growth and development of plant organs 37,38 . The current study showed that endogenous hormones have a leading role from the initial growth stages to the early stages of calli senescence (Fig. 1B). The results from OPLS-DA and the highest value of VIP showed JA to be the dominant endogenous hormone during the initial stages of calli growth, whereas JA and IAA were dominant during the early stages of calli senescence ( Fig. 1B and Table S1). Interestingly, IAA and JA showed similar changes in concentration during the early stages of calli senescence (Fig. 1A), suggesting that they have synergistic effects and jointly inhibit callus growth. The callus with high endogenous IAA content had a higher emergence rate, whereas decreases in IAA content promoted callus dedifferentiation. The content of IAA at 75 d was significantly higher than that at other stages (Fig. 1A), suggesting that the IAA content should be reduced during the later stages of the culture process to enable callus dedifferentiation.
Transcriptome sequencing was performed on different development stages of callus of O. fragrans. In total, 47,340 genes were obtained by transcriptome sequencing, and 1798 previously unidentified genes were found by comparison with the whole O. fragrans genome. These latter genes might be related to the proliferation and differentiation of O. fragrans calli. The KEGG pathway analysis of DEGs showed that these genes were significantly enriched in pathways, including plant hormone signal transduction (Figs. S3, S4). The transcriptome sequencing of maize embryogenic calli showed that DEGs were mainly related to photosynthesis and plant hormone signal transduction 39 . Similarly, 1418 DEGs were identified by transcriptome sequencing from embryogenic and nonembryogenic calli of Picea spruce, mainly involved in plant hormone signal transduction and other pathways 40 . From plant hormone signal transduction pathways, many genes, including some candidate genes that regulate callus proliferation in O. fragrans, were identified. These results provide an interesting direction for the future study of O. fragrans calli.
WGCNA describes the relationship between genes and samples, which can further identify key genes in signaling pathways. Hub genes related to anther development in cotton were identified by WGCNA 41 . Similarly, WGCNA on endogenous hormones and the transcriptome of O. fragrans callus identified modules corresponding to the upregulation of JA and IAA (Fig. 3), revealing five hub genes in two characteristic modules: two hub genes (OfSRC2 and OfPP2CD5) in the JA pathway and three hub genes (OfARR1, OfPYL3, and OfEIL3b) in the IAA pathway (Fig. 5A,B). Endogenous hormones form a complex signaling network in plants, influencing the growth and development of cells by regulating key genes involved in cell proliferation through endogenous signal transmission [42][43][44] . Plants can respond quickly when exposed to external stimuli, resulting in changes in the transcription level of numerous genes. Therefore, a large number of TFs are involved in responding to external stress. Furthermore, we identified that the expression trend of hub genes and TFs by qRT-PCR (Fig. 6A). Finally, we screened out five pairs of combinations with interaction relationships: OfMYC2a and OfPP2CD5, OfERF4 www.nature.com/scientificreports/ www.nature.com/scientificreports/ and OfSRC2, OfHSFA1 and OfEIL3b, OfTGA21c and OfARR1, and OfHSFA1 and OfPYL3 (Fig. 6B). These results imply that hub genes are mainly regulated by TFs, changing the content of JA and IAA, and regulating the proliferation of O. fragrans calli.
In this study, a total of four candidate TFs (OfMYC2a, OfERF4, OfTGA21c, and OfHSFA1) which related to the calli formation were identified and used for the qRT-PCR analysis. However, the homologous genes of the reported calli formation related TFs (bZIP, WOX, and MYB) were not found in the feature modules, which indicated that the gene regulatory mechanism of calli formation in O. fragrans could be different with the model species. MYC TFs are activators or suppressors of JA gene expression, which regulates plant growth and development 45 . Cloning of the AtMYC2 homologous gene MdMYC2 showed that MdMYC2 was involved in the regulation of JA signaling, resulting in the production of more anthocyanins in overexpressed apple calli 46 . ERF TFs occur in plants and are involved in plant responses to biological and abiotic stresses 47 . In Taxus, taxol biosynthesis is regulated by ERF TFs dependent on JA signal transduction 48 . In the current study, early callus growth was rapid (Fig. 7), probably because of TFs OfMYC2a and OfERF4 acting as suppressor factors in response to gene expression, downregulating the expression of hub genes (OfPP2CD5 and OfSRC2) and resulting in decreased JA content and rapid callus growth. TGA TFs are members of the bZIP family and have an important role in the development of stress tolerance 49 . In A. thaliana, TGA TFs enhance resistance by regulating the auxin signaling pathway 50 . HSFA is a heat shock TF involved in regulating the expression of cold, heat, high salt, and other stress-related genes to improve plant stress tolerance 51 . In rice, OsHSFA3 is not only involved in the regulation of plant heat resistance, but also improves plant resistance by regulating the biosynthesis of abscisic acid and reactive oxygen species 52 . The current study also showed that the callus grew slowly or even stopped growing during the early stages of calli senescence (Fig. 7), probably because TFs OfERF4 and OfMYC2a activated the expression of hub genes (OfSRC2 and OfPP2CD5), leading to increases in JA content. At the same time, OfTGA21c and OfHSFA1 also activated the expression of hub genes (OfARR1 and OfPYL3, OfEIL3b), leading to increases in IAA content, inhibiting the growth of O. fragrans calli.

Materials and methods
Plant materials. The plant material for this experiment was O. fragrans 'Rixiang Gui' and the collection of the plant material complied with relevant institutional, national and international guidelines and legislation. Young leaves (the first pair of leaves) of cutting seedlings of O. fragrans 'Rixiang Gui' , which were 2 years old and cultivated in greenhouses at Nanjing Forestry University, were first washed with dishwashing liquid to remove attachments on the surface. They were then washed with pure water in a beaker. The leaves were subjected to  www.nature.com/scientificreports/ dark treatment at 4 °C for 2 h and washed with running water for 1 h. Then, the leaves were transferred to an ultra-clean workbench for disinfection. The explants were subjected to 75% ethanol treatment for 30 s and were rinsed in sterile water three times. Subsequently, they were placed in 5% NaClO for 6 min and then rinsed with sterile water five to six times. The leaf edges and veins were then removed with a sterile scalpel, and the remaining leaves were cut into 0.5 × 0.5 cm pieces and inoculated in induction medium comprising MS (Murashige and Skoog medium, Duchefa, The Netherlands) + 6-BA (6-benzylaminopurine, Biosharp, China) 2.0 mg/L + NAA (naphthaleneacetic acid, Biosharp, China) 0.1 mg/L. The medium was changed every 30 days. Calli of O. fragrans 'Rixiang Gui' at 25 days, 45 days (20 d after the callus began to appear), 55 days (30 days after the callus began to appear), 65 days (40 days after the callus began to appear), and 75 days (50 days after the callus began to appear) from the initial growth stages to the early stages of senescence were removed (Fig. 8) and cut with a sterile scalpel on an ultra-clean workbench. Each callus was transferred to a sterile enzyme-free centrifuge tube with sterile forceps, and immediately frozen in liquid nitrogen for 5 min. All samples were stored at −80 °C. Three samples were collected for each growth period.
Endogenous hormone content measurements. A total of 15 samples across the five growth stages were sent to China Agricultural University for measurement of endogenous hormone contents. The callus (0.5 g) of each treatment were well-powdered using liquid nitrogen and then samples were crushed by cold methanol. The extract was achieved using 2 mL of 80% methanol at 4 °C for 4 h. Then, the samples were centrifuged for 10 min at 6000 rpm/min to absorb the supernatant. The supernatant was extracted on a C-18 SPE column and evaporated to dryness in a rotary evaporator. After that, the samples were resuspended in 1 ml of sample diluent RNA extraction, cDNA library preparation, and sequencing. Total RNA was extracted from the 15 samples using an RNA Purification Kit (Invitrogen, Carlsbad, CA, USA). The integrity of the RNA was verified by RNase-free agarose gel electrophoresis and the concentration was measured using a Nano Drop 2000 spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA). The enriched mRNA was then fragmented into short fragments using fragmentation buffer and reverse transcribed into cDNA with random primers. Second-strand cDNAs were synthesized by DNA polymerase I, RNase H, dNTP, and buffer. The the cDNA fragments were then purified with a QiaQuick PCR extraction kit (Qiagen, Venlo, The Netherlands) and end repaired, poly(A) was added, and the fragments were ligated to Illumina sequencing adapters. The ligation reaction was purified with the AMPure XP Beads (1.0X). Ligated fragments were subjected to size selection by agarose gel electrophoresis and amplified by polymerase chain reaction (PCR). For each stage of O. fragrans callus development, three RNA samples were used to construct a cDNA library and for Illumina Novaseq6000 sequencing, which was completed by Gene Denovo Biotechnology Co. (Guangzhou, China). All sequenced genes were blasted with the genome of O. fragrans.
De novo assembly of RNA-Seq reads and quantification of gene expression. First, the adapter sequences were deleted from the raw reads. Then, low-quality reads (with over 40% of bases with quality scores of 10 or lower and/or over 10% of bases unknown) were deleted by using fastp version 0.18.0 (https:// github. com/ OpenG ene/ fastp) from each data set to establish more reliable results. Following this, clean, high-quality reads from all the samples were combined. Then an index of the reference genome was built, and paired-end clean reads were mapped to the reference genome using HISAT 2.2.4 (https:// daehw ankim lab. github. io/ hisat2/) with '-rna-strandness RF' and other parameters set as a default. The mapped reads of each sample were assembled by using String Tie version 1.3.1 (http:// ccb. jhu. edu/ softw are/ strin gtie/) in a reference-based approach 54,55 . For each gene, the expression level was measured by fragments per kilobase exon model per million mapped fragments (FPKM), based on the number of uniquely mapped reads, to eliminate the influence of different gene lengths and sequencing discrepancies in the gene expression calculation 3 .
Identification and functional analysis of DEGs. Transcripts with a fold change > 2 and false discovery rate (FDR) < 0.05 were considered to be differentially expressed between growth stages. The gene sequences were compared with GO and KEGG databases [56][57][58] . The threshold was set to E-value < 1e−5, and the gene function information was obtained 59,60 . GO enrichment analysis and KEGG pathway enrichment analysis, the former being a standardized gene function classification system and the latter the main pathway analysis database, were performed for DEGs 61 . The annotations identified the signal transfusion pathways and related metabolic pathways of DEGs, which further revealed the biological functions of these genes 62 .
Screening of hub genes and co-expression network regulation graphs. WGCNA was used to analyze the correlation between endogenous hormones and DEGs. The endogenous hormone content in O. fragrans was measured as a physiological indicator. A gene clustering tree was constructed according to the correlation between the expression levels of genes and the gene module. Modules were divided according to the clustering relationship between genes. The modules with similar expression patterns were combined according to the similarity of the module feature values. Genes with similar expression patterns were placed in the same www.nature.com/scientificreports/ module. We first determined the module most related to the character of interest. Based on KEGG annotation of the transcriptome, combined with the correlation analysis diagram of traits, corresponding feature modules were found. The MM (MM > 0.9) value of the characteristic module corresponding to key genes was set to further identify hub genes. To study the relationship between hub genes and edge genes, Cytoscape software was used to draw a regulatory network map. First, the downstream genes were sorted according to the weight value from the highest to the lowest. Then, according to the weight value (weight > 0.2), the TFs were selected to draw the regulatory network map.
Verification of gene expression using qRT-PCR. To verify the accuracy of the transcriptome data, genes were selected from transcriptome data. Primers were designed using Primer 5.0 software (Table S3).
OfRAN of O. fragrans was selected as the internal reference gene 63 . cDNAs were synthesized from 5 µg total RNA and diluted 20-fold for gene expression experiments. The qRT-PCR reaction kit used was SYBR Premix Ex Taq (Takara Biotechnology, Dalian, Liaoning Province, China). Reactions were performed at 95 ℃ for 3 min, followed by 40 cycles of 95 °C for 5 s, and 60 °C for 30 s. Relative gene expression levels were calculated according to the 2 −∆∆Ct comparative CT method 3 . Three technical replicates and three biological replicates were used for each sample.

Conclusion
This is the first study of the transcriptome and physiology of O. fragrans callus. The content of endogenous hormones and transcriptome sequencing of O. fragrans calli were studied at different growth stages. The initial stages of callus growth were mainly regulated by JA, whereas the early stages of callus senescence were regulated by JA and IAA. The association between the co-expressed gene modules formed by DEGs and the endogenous hormones was analyzed by WGCNA. Here, we found that hub genes (OfSRC2, OfPP2CD5, OfARR1, and OfPYL3, OfEIL3b) involved in plant hormone signal transduction pathways could be regulated by four TFs (OfERF4, OfMYC2a, OfTGA21c, and OfHSFA1), which altered the JA and IAA content and regulated callus proliferation in O. fragrans. These results provide new insights for future research on the proliferation efficiency of O. fragrans callus.

Data availability
The transcriptome data have been uploaded to the NCBI Sequence Read Archive (https:// www. ncbi. nlm. nih. gov/ sra/) under accession number PRJNA760274. The datasets used and/or analyzed during the current study are available from the corresponding author on reasonable request.